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1 Introduction 

Most massive stars (8 M Q — » 80 M ) must transition at the ends of their lives 
into neutron stars or stellar-mass black holes. That they do so when their low- 
entropy cores reach the Chandrasekhar mass, gravitationally collapse, and launch 
a supernova explosion has been demonstrated both by direct observations (cf., 
the neutrinos from SN1987A) and by a host of compelling theoretical arguments. 
However, numerical simulations of the process of core collapse, bounce near 
nuclear densities, shock wave generation, and shock propagation have failed to 
recreate in detail the observed gravitational masses of known neutron stars, the 
expected nucleosynthesis pattern, and empirical supernova energies. All semi- 
realistic, one-dimensional (spherically-symmetric) simulations conducted to date 
hzzle into quasi-static accreting proto-black holes. What is more, they seem to 
do so convincingly, despite concerted attempts over the years to include all the 
known neutrino and nuclear physics or general relativity |]^,poP,|5l, |l^ , |l^|Jr7|^] . 
It is now fairly clear that the devil is not "in the details" and that ID models 
don't explode. Something large and major, not something at the "10-20% level," 
would have to be missing to alter this conclusion. 

Two-dimensional simulations conducted in the nineties [pr|J^ Jl5| , p^| and more 
recent three-dimensional SPH simulations [fl3| do explode (when the correspond- 
ing simulations in one dimension do not), but these multi-D numerical experi- 
ments do not produce supernova explosions that satisfy all of the above obser- 
vational constraints. Furthermore, they all employ some realization of a simple 
flux-limited, energy-integrated diffusion algorithm, in lieu of full neutrino trans- 
port (which is quite difficult in 2- or 3-D). It is thought by many that the 
apparent marginality of success demands that the transport be handled with 
a multi-group, multi-angle technique (and that the hydrodynamics be handled 
in full general relativity) before the multi-D explosion simulations are to be 
considered valid, or even indicative of the mechanism. Perhaps. However, this 
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puts too much reliance on hardware and software, and not enough on physical 
understanding. 

So, what is wrong? And how to fix it? In order to see where we are going we 
must first see where we have been. 



2 One-Dimensional Simulations: The Spherically- 
Symmetric Paradigm 

The primary elements of the neutrino-driven mechanism were mapped out in 
the work of Colgate and White |Tl[] and ArnettQ. In the former, a cold (low- 
entropy), already deleptonized core was collapsed and after bounce at nuclear 
densities a neutrino luminosity which heated the outer mantle was turned on for 
14 milliseconds. The neutrino transport was artificial, but it drove an explosion 
with 5 x 10 52 ergs. In the latter, a higher-entropy core was collapsed, which 
bounced due to thermal nucleon pressure at a central density near a tenth of 
nuclear, and neutrino transport was handled more realistically with an energy- 
integrated diffusion approach. Neutrino diffusion was found to transport energy 
from the hot interior to the mantle on close to the dynamical timescale of the 
newly-generated shock wave. In Arnett's calculations, this kept the shock wave 
from stalling, an explosion ensued, and the total duration of the neutrino lumi- 
nosity was less than 100 milliseconds, but longer than in Colgate and White. 
In both cases, neutrinos played the central role in launching the explosion and 
a neutron star remained, but the initial progenitor models, neutrino transport, 
neutrino physics, and equations of state were primitive by today's standards. 
The timescale for the cooling phase (detected for SN1987A) was off by a factor 
of 100-1000, the spectra were much harder than current simulations indicate, and 
neutrinos were not included. Nevertheless, the core idea of the neutrino- 
driven mechanism, that neutrinos generated in the inner core and deposited in 
the outer mantle drive the explosion, was implicit in both papers. 

The essential ingredient of the neutrino-driven mechanism is the transfer of 
gravitational binding energy from the core to the mantle of the newly-formed 
protoneutron star through the mediation of neutrinos liberated at the high tem- 
peratures and densities generated in its interior. An efficient coupling is required, 
that, if not achieved, will lead to a fizzle. With up-to-date nuclear equations of 
state, neutrino physics, progenitor models, and transport algorithms, the cou- 
pling efficiency in one-dimensional simulations has been determined to be inad- 
equate both to forestall the stagnation of the bounce shock and to reignite it 
once it has stalled. Figure [l] depicts radius trajectories of mass zones versus time 



for one of our ID multi-group, Feautrier/tangent-ray calculations (§2.2). The 
progress of the shock and its stagnation during the first 250 milliseconds after 
bounce are manifest. Wilson [^6| and Bethe and Wilson^) suggested that shock 
stagnation is only temporary, leading to the "delayed" scenario, but without 
the neutron-finger convective boost (an inherently multi-D effect which others 
have failed to reproduce) in the driving neutrino luminosity Wilson himself did 



Neutrino-Driven Supernova Mechanism 3 



not find explosions. In sum, embellishments with general relativity, the employ- 
ment of multi-group, multi-angle Boltzmann solvers, and refinements in neutrino- 
matter interactions have not changed the conclusion that in spherical symmetry 
the coupling efficiency of the emergent neutrinos to the protoneutron star man- 
tle is too small to lead to explosion, even after a pause P,p5P,p|, ^8| , p^l7|p^ | . 
Unfortunately, this conclusion seems to be independent of progenitor model. 

2.1 Supernova Energetics Made Simple 

It is important to note that one is not obliged to unbind the inner core (~10 
kilometers) as well; the explosion is a phenomenon of the outer mantle at ten 
times the radius (50-200 kilometers). One consequence of this goes to the heart of 
a general confusion concerning supernova physics. Though the binding energy of 
a cold neutron star is ~3 x 10 53 ergs and the supernova explosion energy is near 
10 51 ergs, a comparison of these two numbers and the large ratio that results 
are not very relevant. More germane are the binding energy of the mantle (inte- 
rior to the shock or, perhaps, exterior to the neutrinospheres) and the neutrino 
energy radiated during the delayed phase. These are both at most a fewxlO 52 
ergs, not ~3 x 10 53 ergs, and the relevant ratio that illuminates the neutrino- 
driven supernova phenomenon is ^10 51 ergs divided by a fewxlO 52 ergs. This 
is ~5-10%, not the oft-quoted 1%, a number which tends to overemphasize the 
sensitivity of the neutrino mechanism to neutrino and numerical details. 



Fig. 1. Radius trajectories of selected mass zones in our 300-zone, 20 energy-group 
calculation of stellar collapse, bounce, and shock propagation for the Woosley and 
Weaver 11 Mq progenitor^]. In this simulation, as in all other ID simulations, the 
shock does not revive. 
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Furthermore, there is general confusion concerning what determines the su- 
pernova explosion energy. While a detailed understanding of the supernova mech- 
anism is required to answer this question, one can still proffer a few observations. 
First is the simple discussion above. Five to ten percent of the neutrino energy 
coursing through the semi-transparent region is required, not one percent. Im- 
portantly, the optical depth to neutrino absorption in the gain region is of order 
~0.1. The product of the sum of the v e and v e neutrino energy emissions in the 
first 100's of milliseconds and this optical depth gives a number near 10 51 ergs. 
Furthermore, the binding energy of the progenitor mantle exterior to the iron 
core is of order a fewxlO 50 to a fewxlO 51 ergs and it is very approximately this 
binding energy, not that of a cold neutron star, that is relevant in setting the scale 
of the core-collapse supernova explosion energy. Given the power-law nature of 
the progenitor envelope structure, it is clear that this binding energy is related 
to the binding energy of the pre-collapse iron core (note that they both have 
a boundary given by the same GM/R), which at collapse is that of the Chan- 
drasckhar core. The binding energy of the Chandrasckhar core is easily shown 
to be zero, modulo the rest mass of the electron times the number of baryons in 
a ~1.4 M Q Chandrasckhar mass. (The Chandrasekhar mass/instability is tied 
to the onset of relativity for the electrons, itself contingent upon the electron 
rest mass). The result is ~10 51 ergs. 

The core-collapse explosion energy is near the explosion energy for a Type 
la supernovae because in a thermonuclear explosion the total energy yield is 
approximately the 0.5 MeV/baryon derived from carbon/oxygen burning to iron 
times the number of baryons burned in the explosion. The latter is >half the 
number of baryons in a Chandrasekhar mass. The result is ~10 51 ergs. This 
is the same number as for core-collapse supernovae because 1) in both cases 
we are dealing with the Chandrasekhar mass (corrected for electron captures, 
entropy, general relativity, and Coulomb effects) and 2) the electron mass and 
the per-baryon thermonuclear yield are each about 0.5 MeV. 

While more detailed calculations are clearly necessary to do this correctly, the 
essential elements of supernova energetics are not terribly esoteric (if neutrino- 
driven), at least to within a factor of 5, and should not be viewed as such. 

2.2 Some Recent Results using our Feautrier/Tangent-Ray/ALI 
Method 

One is firm in the general conclusion that ID models with good physics and 
numerics do not explode because, coming from many different angles, various 
groups have now verified this. We here briefly describe our contribution to this 
activity. 

It is thought that the angular distribution of neutrinos in the gain region 
needs to be calculated with precision, since neutrino energy deposition is pro- 
portional to neutrino energy density and, for given luminosities/fluxes, this is 
sensitive to v e and v e angular distributions that can not be derived using flux- 
limited diffusion. Furthermore, the stiff dependence of the absorption cross sec- 
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tions on neutrino energy requires multi-group approaches. Hence, multi-group, 
multi-angle Boltzmann transport algorithms are to be preferred. 

We have constructed such a code, using the Feautrier variables advocated in 
standard stellar atmospheres work|l9) and the tangent-ray method to establish 
a dense angular gridQ. Our transport solver calculates in the comoving frame, 
uses accelerated A iteration ( ALI) j|Jl^,|2^,^l| to speed convergence of the solu- 
tion, is implicit in time, second-order accurate in space, and iterates between 
the Boltzmann/transfer equation and the zeroth- and first-moment equations 
(in a method akin to the variable Eddington factor approach) until a converged 
global solution to the full transport equation is achieved. No ad hoc flux limiters 
or artificial closures are necessary. This iteration scheme is fast (convergence to a 
part in 10 6 in 2 to 10 steps) and automatically conserves energy in the transport 
sector. The Feautrier scheme can transition to the diffusion limit seamlessly and 
accurately and the tangent-ray method automatically adapts with the hydrody- 
namic grid as it moves. In constructing the tangent rays, we cast them from every 
outer zone to every inner zone. Hence, if there are 200 radial zones, the outer 
zone has 199 angular groups. Because of the spherical nature of the core-collapse 
problem and the need to accurately reproduce the angular distribution of the ra- 
diation field that transitions from the opaque (inner) to the transparent (outer) 
regions, such fine angular resolution is useful, though computationally demand- 
ing, as radiation becomes more and more forward-peaked. Figure depicts a 
snapshot of the angular distribution of the f e neutrinos at various energies from 
1 to 320 MeV, at a radius of 42 kilometers, 40 milliseconds after bounce. The 
positions of the angular bins are shown as dots in the lower hemisphere. The 
progressively more forward-peaked distribution at lower energies (less coupled) 
is clearly well-resolved. The code calculates the Feautrier variables (and, hence, 
the specific intensity) at every radial zone, for every energy group, for each neu- 
trino species (we follow the standard "3"), at each timestep. We employ 20-40 
energy groups either logarithmically or linearly spaced from 1 MeV to either 100 
MeV (for v e and V M ") or 320 MeV (for v e ). 

The hydrodynamics is explicit in time and Lagrangean, uses a predictor/corrector 
method, employs artificial viscosity to handle shocks, and is Newtonian. We had 
originally used a PPM hydro scheme, but could not easily incorporate neutrino 
radiation pressure into the Riemann solver. We have constructed a dense equa- 
tion of state table in the variables T, p, and Y e that contains the 14 variables 
needed by the radiation hydrodynamics. We have constructed additional tables 
of the z/j-electron scattering kernels and the e + /e~ annihilation kernels d^5|. 
Neutrino-electron redistribtion/scattering is handled explicitly, which is a great 
time and memory saver. The radiation/matter couplings (both the energy and 
electron fraction updates) are handled implicitly and in operator-split fashion, 
using the ALI to facilitate convergence. Figure || provides snapshots we obtain of 
a representative evolution of the Y e profile before and after trapping and bounce. 

The minuses of our approach are that it is Newtonian, that the hydro is 
explicit, and that, being Lagrangean and not adaptive, the envelope must be 
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pre-zoned densely in mass to maintain reasonable resolution of the shock at late 
times. Rampp and Janka]2^,^3| use an Eulerian grid and PPM hydro, remap be- 
tween the comoving frame radiation solution and their static hydro grid, employ 
a variable Eddington factor/tangent-ray/Feautrier scheme as well, and operator 
split the neutrino-matter couplings. They cast tangent rays to a subset of interior 
Eulerian zones and thereby maintain a static, fixed set of angles. In addition, 
they have a simple method for approximately incorporating general relativity. 
There are lots of minor differences between our numerical implementations, but 
despite them our simulation results are quite similar fl25f. 



Liebendorfer et al. |L6 17 have assembled a very different radiation/hydrodynamic 
code complex that has both strengths and weaknesses. Among their many strengths 
is that their code is fully implicit, adaptive (a great advantage), and general rel- 
ativistic. Minor weaknesses are that they operator-split all the terms in the 
transport equations separately and use perforce (due to their implicit redistri- 
bution method) a small number of angles (6) in the S„ method they utilize to 
handle the angular distribution of the radiation field. S n can not handle forward- 
peaked radiation fields, but for the core collapse problem, the gain region inte- 
rior to a stalled shock (where the angular distribution is most problematic) is 
sufficiently compact that their S„ technique is more than adequate. Where it 
fails, in the outer regions, the neutrino/matter coupling is not so germane to 
the problem of the supernova mechanism. There are many other differences of 
implementation and approach, but, again, despite these differences their results 



Fig. 2. A snapshot of the angular distribution of the v e specific intensity at a radius 
of 42 kilometers, 40 milliseconds after bounce, for energy groups from 1 to 320 MeV. 
The 11 Mq progenitor of Woosley and Weaver |27j was used. 
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for velocities, trapped electron fractions, entropies, densities, and the effect of 
neutrino-electron redistribution are both qualitatively and quantitatively similar 
to those of both ourselves [^5) and Rampp and Janka[ p2[|2^ |. 

Figure ^ portrays the evolution of the luminosity of the electron neutrinos 
with time up to 200 milliseconds after bounce for three progenitor models. The 
breakout peak, and the precursor peak just a few milliseconds earlier, are not 
only similar to one another, but they are similar to the corresponding results 
of Rampp and Janka^2|. Though mostly similar to the results of Liebendorfer 
et al. off the peak, their peak luminosities can be 20-30% higher than ours. 
Due to the 20-M Q model's thicker envelope, the subsequent accretion luminosity 
for this progenitor is ^60% higher than that for the two other progenitors 200 
milliseconds after bounce. In none of our ID simulations, carried out to as long as 
~1 second after bounce and with the best neutrino and nuclear physics available, 
do the models explode. Moreover, we have incorporated, in approximate fashion, 
the possible effects of the neutrino oscillations, using the numbers derived from 
the solar and atmospheric neutrino experiments, and have seen no effect. This is 
not unexpected, since the "matter effect" severely suppresses oscillation between 
flavor neutrinos. 

We are forced to conclude, on the basis of both our simulations and those 
of others, that the program to determine whether core-collapse supernovae ex- 
plode in spherical symmetry (ID) when all the best physics and numerics are 
implemented has failed by a comfortable (?) margin. This leads us of necessity 
to multi-D effects. 



Fig. 3. A collection of snapshots of the Y e profile versus interior mass for the Newtonian 
evolution of the Woosley and Weaver p7| 11 Mq progenitor. A total of ~250 millisec- 
onds is depicted and the simulation was carried to ~50 milliseconds after bounce. 



8 



Adam Burrows and Todd A. Thompson 



3 Multi-Dimensional Simulations 

The success of the 2D BHF[Q and Fryer et al.|l2| simulations can be traced 
to the increase in the efficiency of the neutrino-matter coupling in the gain re- 
gion interior to the shock due to the longer dwell time of convecting/absorbing 
parcels of matter. In ID, accreted matter falls straight through the gain region. 
Consequently, in ID the increase in entropy due to neutrino energy deposition 
is modest. In 2- or 3-D, matter resides a bit longer in the gain region due to 
convection. As a result, the steady-state entropy of an average parcel of matter 
is higher. This translates into a larger gain region, that is on average less bound 
and more unstable to favorable changes in the global conditions of accretion and 
luminosity. In the language of Burrows and Goshy||, the global instability con- 
dition of the protoneutron star is eased and changes in the control parameters 
(such as M and ) can more easily result in explosion. Whether they do when 
better neutrino transfer and general relativity are implemented is unclear, but 
recently Janka and Rampp (JR, this volume) have seen a very weak explosion 
in a calculation using 2D hydro in a 27° wedge and their ID Feautrier transport 
along rays. They obtain no explosion in ID for the same inputs. In this 2D cal- 
culation, an average hydro solution in the wedge is used as input to a spherical 
calculation of the Eddington factors, and these are employed for all the individ- 
ual Feautrier solutions along rays. This is better than the flux- limited solution 
along rays of BHF, incorporates better neutrino transport, has an approximate 
prescription for relativity, but seems to result in a weaker explosion. However, 
the JR algorithm is not yet true 2D transport and the calculations need to be 



Fig. 4. The electron-neutrino luminosity versus time around bounce for three different 
progenitors. The primary focus of this plot is the prodigious shock breakout burst. 
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done over a full 180° angular region. Nevertheless, it is clear from the JR calcu- 
lation and those of BHF[Q], Herant et al.[[l4], and Fryer et al.|L2| that convection 
makes the protoneutron star mantle more unstable. What is not clear is whether 
this is enough. 

Fryer and Warren) 13 have recently demonstrated that their 3D SPH simula- 
tions explode quantitatively and qualitatively in much the same way as their 2D 
simulations, despite a general concern that 2D and 3D convection have differ- 
ent clump size spectra and cascade character. This seems not to matter much; 
it is the overall global behavior of the convecting region that matters and, ac- 
cording to them, the largest size and velocity scales. Small scales are a detail. 
The former are set by the size of the unstable region and are similar in both 
2D and 3D. However, these calculations need to be verified with multi-group, 
multi-angle techniques, first in 2D (which has yet to be done), then in 3D (a 
major computational challenge for the future). 



4 Coda 

Multi-dimensional effects seem to be required to ensure or enable core-collapse 
supernova explosions. However, even though what the essential elements for ex- 
plosion are is unclear, we believe that the solution, when found, will not be 
marginal. A number of groups are now embarked upon separate development 
paths to credible 2D (later 3D) radiation/hydro schemes (Arizona/Israel, ORNL, 
MPA). We, in collaboration with Eli Livne and Itamar Lichtenstadt, are devel- 
oping a 2D, moving grid, discontinuous finite clement (DFE), multi-group, and 
multi-angle (S„) neutrino radiation/hydrodynamics code. With it, we hope soon 
to simulate multi-D collapse and to test out the numerous avenues now before 
us. Stay tuned. 
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